Calibration of camera internal parameters based on grey wolf optimization improved by levy flight and mutation

Traditional calibration technology has been widely used in measurement and monitoring; however, there are limitations of poor calibration accuracy, which can not meet the accuracy requirements in some scenarios. About this problem, we proposed a grey wolf optimization algorithm based on levy flight and mutation mechanism to solve camera internal parameters in this paper. The algorithm is based on the actual nonlinear model, which takes the minimum average value of reprojection error as the objective function. The grey wolf position is randomly generated within a given range. Then, the grey wolf optimization algorithm based on levy flight and mutation mechanism is used to iteratively calculate the optimal position, which is the internal parameters of cameras. The two groups of experimental data were performed to verify the algorithm. The result shows better effectiveness and calibration accuracy of the proposed algorithm compared with other optimization methods.


Scientific Reports
| (2022) 12:7828 | https://doi.org/10.1038/s41598-022-11622-w www.nature.com/scientificreports/ those methods to camera calibration, which can obtain the accuracy of camera parameters. Li, et al. 8 proposed an effective camera calibrate method, which combined differential evolution way with PSO algorithm based on the camera calibration model of geometric parameters and lens distortion. Li, et al. 9 proposed to combine the genetic algorithm and PSO algorithm to solve camera calibration problem, which can avoid particle swarm to the local optimal. Liu, et al. 10 proposed an improved PSO algorithm with a new evaluation function, which adjusted the inertia factor and learning factor. This method achieved good results in camera internal parameter calibration. Qin et al. 11 proposed a full-parameter PSO algorithm based on mutation strategy, improving convergence speed and accuracy compared with the traditional particle swarm calibration ways. Yang, et al. 12 utilized the simulated annealing algorithm and the PSO algorithm to combine the advantages of those methods through a collaborative mechanism, which the hybrid algorithm can merge the advantages of different methods. the method produces more good performance than a single algorithm, which can improve the accuracy of the calibration results. The preceding observation 13 suggests that the results of calibration is affected by the accuracy of feature extraction. Thus, numerous researchers have concentered on creating targets with distinctive features that can be accurately localized in the images, including 1D, 2D and 3D targets. Among those targets, 3D targets are difficultly created and hardly solved. Therefore, it limits their applications. The advantage of 1D targets is flexibility, however, due to the fewer feature point involved, they are hard to acquire high accuracy. In contrast, the second objectives have been broadly investigated attributed to their flexibility and effectiveness. There are two typically used patterns of 2D targets: checkerboards and circles. However, in some scenarios, the calibration accuracy of 2D targets is still limited. In this paper, an improved grey wolf algorithm based on Levy flight and mutation strategy is designed to solve camera calibration; and then the objective model is optimized by this method to acquire the accuracy parameters of camera. The error of calibration parameter's results obtained is smaller than the grey wolf algorithm, particle swarm optimization algorithm, and Zhang's calibration method 7 .

Monocular camera model and four coordinate systems
Camera imaging model, which mainly reflects the process of the camera taking pictures of the actual 3D world. The model of camera imaging is the foundation of camera calibration and the solution camera parameters. According to the difference of model of camera imaging, it can be categorized into linear and nonlinear model. The linear model is based on the pinhole imaging precept, which establishes the relationship of pinhole imaging point and the corresponding object point through its geometric courting 14 . The camera distortions are usually caused by the physical structure in practical applications. Therefore, a nonlinear model may be applied to correct the camera distortion. The cameral imaging model is shown in Fig. 1.
Imaging includes four coordinate systems. They are called World Coordinate System (WCS), Camera Coordinate System (CCS), Image Coordinate System (ICS) and Pixel Coordinate System (PCS). They are converted according to a certain scale.
World coordinate system (X W , Y W , Z W ) , which is also called measurement coordinate system. It is an orthogonal 3D rectangular coordinate system, which is established based on an object in reality.(X WP , Y WP , Z WP ) denotes the coordinate of object point P in the world coordinate system.
Camera coordinate system (X C , Y C , Z C ) : It is a 3D rectangular coordinate system. The origin is located in the optical center of the camera. The axis Z C is the main optical axis of the camera.(X CP , Y CP , Z CP ) is the image coordinate of object point P in the ideal pinhole model.
Image coordinate system x, y : A 2D rectangular coordinate system established on the image plane. The transformation between the CCS and ICS is known as perspective transformation. The coordinate axis x and y  www.nature.com/scientificreports/ are respectively parallel to the axis X c and Y c .The origin O 1 is situated in the intersection of the camera optical axis and the CCD image plane. x P , y P denotes the coordinate of object point P in the image coordinate system. Pixel coordinate system (u, v):It is established on the same plane as the image coordinate system, while the difference is that the origin of the two is different. The coordinate axis x and y are respectively parallel to the axis u and v. u P, v P denotes the pixel coordinate of object point P in the pixel coordinate system. External parameters of the camera. Rigid body transformation describes the relative positional relationship between the world coordinate system (WCS) and the camera coordinate system (CCS) through the external parameters of the camera. Assume that there is a point P in three-dimensional space in reality, and its coordinate in the world coordinate system is (X WP , Y WP , Z WP ) . The formula (1) shows the coordinates (X CP , Y CP , Z CP ) of point P when it's converted to camera coordinates.
where R and T respectively represent the rotation matrix and translation vector. And they do not depend on the inherent attributes of the camera itself. It shows there are two external parameters involved in camera calibration.
Internal parameters of the camera. If the error caused by distortion is not taken into account, as shown in Fig. 1, camera model is a linear model. And then the geometric knowledge can be used to solve the relationship between parameters. The point coordinates of camera (X CP , Y CP , Z CP ) is calculated by external parameters. After a perspective transformation, the points on the physical coordinate system of the image are obtained as x P , y P .
where f is the camera's focal length.
Imaging in the image pixel coordinate system is the last presentation form of the camera. Pixel coordinates (u, v) can be obtained directly from a photograph. Also taking the point P as an example. The transformation between the pixel point and the image physical coordinate system in the same imaging plane with a different origin can be described by the following formula: where (u 0 , v 0 ) is the intersection point of the two coordinate axes, namely the coordinates of the origin O 0 , d x and d y represents the physical size of the unit pixel on the two coordinate axes.
In combination with (1), (2) and (3), it can be obtained that: where K and F are the internal parameter matrix and the external parameter matrix respectively.
Solving camera parameters. So as to simplify the calculation, it's miles assumed that the template plane is on Z w = 0 in WCS. Formula (4) can be simplified as below: In order to solve camera internal parameters, a new concept called homography matrix is introduced. Truly understood as it is used to explain the position mapping relationship between the object within the world coordinate system and the pixel coordinate system. The corresponding transformation matrix is the homography matrix, which is expressed as: Then, the expression can be expressed as: Because the image point in the central area of the image is less affected by distortion factors. Therefore, in order to reduce the calculation error, select the image center area n (n > 9) to pixel points (u P , v P ) and (X WP , Y WP , 0) , the matrix H can be calculated. According to the definition of H in formula (6), r1 and r2 are orthogonal vectors, thus: Namely: where B = K −T K −1 ,and is a symmetric matrix, which can be defined as a 6- The camera intrinsic parameters are determined by the two main constraints, which can be expressed as follows: . In Eq. (11), N ij is known; The matrix B contains 6 unknowns and requires 6 equations to solve. One homography matrix H can express two equations, for which three calibration pictures with different shooting angles are needed to solve B . Then, matrix B is decomposed by Cholesky to obtain the camera's internal parameter matrix K . The specific solution method can be referred to paper 8 . According to the homography H and intrinsic parameters, camera extrinsic parameter can be solved as follows: where the scale factor s is determined by the orthogonal condition. So far, the initial internal and external parameters of the camera have all been solved.
Distortion correction. Note that lens distortion of the camera became no longer considered throughout the above calibration system. The distortion resulting from the lens shape is referred to as radial distortion. Within the pinhole model, a straight line projected onto the pixel plane remains a straight line. however, in actual photos, the digicam lens regularly turns a straight line within the actual environment into a curve in the picture. The closer to the edge of the picture, the more obvious this phenomenon. For the reason that actual processed lens is often symmetrical, irregular distortion is generally radially symmetrical. They specifically divided into two categories, barrel distortion and pincushion distortion.
Barrel distortion is because the photo magnification decreases as the distance from the optical axis will increase, while pincushion distortion is simply the opposite. In those sorts of distortions, the instantly line    X WP = u P h 11 +v P h 12 +h 13 u P h 31 +v P h 32 +h 33 Y WP = u P h 21 +v P h 22 +h 23 u P h 31 +v P h 32 +h 33 www.nature.com/scientificreports/ passing through the intersection of the picture center and the optical axis can keep its form. Similarly, to the radial distortion brought via the form of the lens, tangential distortion can also be delivered due to the inability to make the lens and the imaging surface strictly parallel at some point of the meeting system of the digicam 15 .
In summary, distortion is mainly classified into three types: eccentric, radial and thin prism distortion. Under the comprehensive action, there will be errors in both radial and tangential directions, and the error calculation formula is as follows: where k 1 , k 2 , p 1 , p 2 , k 3 are the five variables, which are the factors related to the distortion correction in both directions.
Here the nonlinear model introduced distortion is used to calibrate the camera.
In summary, the distortion factor k 1 , k 2 , p 1 , p 2 , k 3 , the parameters related to the focal length f x , f y , and the image center (u 0 , v 0 ) are the nine internal parameters of the camera in total.

Improved grey wolf optimization algorithm
Wolf optimization algorithm. Grey Wolf Optimization Algorithm (GWO) is a recent population-based bionic algorithm, which is given by Mirjalili 16 . This algorithm simulates social hunting conduct of grey wolves. In comparison with different meta-heuristic algorithms, the GWO algorithm can be definitely understood as an efficient optimization technique that simulates the social behavior and management hierarchy of grey wolves. In the mathematical instance of the hierarchy of the grey wolves, alpha(α) wolf is referred to as the most suitable solution. consequently, beta(β) wolf is referred to as the second one most appropriate solution after which the following appropriate solution is known as delta(δ). the alternative populations which constitute the furthest solutions are seemed to be the omegas (ω) 17 . The quality wolves must be dealt with as α, β, and δ that help exceptional wolves (ω) in exploring more favorable regions of solution area, as shown in Fig. 2.
The grey wolves encircle the victim while hunting. The surrounding strategy, is mathematically represented as in (14): − → X p (k) represent location of wolves and location of the prey respectively.k shows recent iteration.a = 2 − k k max is a temporal parameter that its elements are linearly decreased from 2 to 0. r 1 and r 2 are random vectors between [0, 1].
These animals need to discover and enclose the location of the quarry. The searching technique is usually directed through the alpha types. In a few conditions, the beta and delta would possibly contribute in searching as well. In GWO, it has been intended that the alpha (best solution), beta and delta can guesstimate the in all likelihood scenario of the sufferer. Hence, the primary three best solutions should be recorded to conduct other hunters. It's miles clear that the end position is located in a random vicinity inside a hypersphere within the search space. Therefore, elite solutions will approximate the top of the line and other hunt retailers simply ought to revise their locations around a predicted victim based totally on stochastically. www.nature.com/scientificreports/ After encircling the prey, the wolves hunt the prey under the leadership of the alphas, betas and deltas 18 . The position of the prey also changes as it runs away. At this time, the positions of the wolves change accordingly, and the calculation formula is as follows: Levy flight. Levy flight is a type of walk that follows Levy probability distribution function and arbitrarily orientations the step length 19 . In addition, Levy flight has been found to be widespread in nature. Creatures such as sharks and cuckoos move according to this law. This distribution is a simple power -law formula. The formula is usually as follows: where s is step length, µ is parameter is location or shift parameter,,γ , L(s) represents the distribution of s.γ > 0 parameter is scale (controls the scale of distribution) parameter. Levy distribution can be represented by a clear power-law equation as: In the Mantegna strategy, the normal distribution is used to solve the random step size, and the method is as follows: where u and v are drawn from normal distributions.
The result of levy flight simulation is shown in Fig. 3. It can be obviously seen from the figure that this random walk can produce high-frequency short step lengths and intermittent long step lengths. This feature can suppress the shortcomings of the grey wolf algorithm that is not difficult to fall into local extremes, and will not weaken the ability of global optimization.

Grey wolf optimization improved by Levy flight and mutation mechanism.
When using swarm intelligence algorithms to solve complex high-dimensional global optimization problems, balancing algorithm convergence and population diversity is very important. As a swarm intelligence algorithm, GWO also is required to solve this problem. Mutation can promote the evolution of living organisms and bring next-generation individuals more adaptable to the environment. The bionic structure of the GWO algorithm gives it the advantages of realization and flexibility. However, the diversity of individuals in the later stage of GWO is insufficient, and the range of prey hunting will be stagnant in the local Optimal area 20 . To maintain the diversity of the following individuals and jump out of the local optimal area, this paper generates grey wolf intermediates according to the mutation mechanism in the differential evolution algorithm, as showed in Eq. (21): where G i (k + 1) is intermediate,i = r 1 = r 2 , X α is the optimal solution, X r2 and X r3 are 2 position vectors randomly selected in a population.F is the scaling factor, 0.5 here.
In the grey wolf optimization algorithm, the state of the grey wolf can be judged by A : when A > 1 , the wolves are far away from the prey, global search is conducted; when A < 1 , the wolves are close to the prey, www.nature.com/scientificreports/ local search is conducted. Therefore, this paper is based on the formula (13), in order to balance the local and Global optimization capability, combined with Levy flight, the formula is as follows: where ε is the step control factor, here is 0.01. The intermediate is obtained through the above formula and then compared with the original individual, and the optimal position is updated according to the greedy (GS) algorithm. www.nature.com/scientificreports/ The optimal individual is chosen at each iteration to achieve optimization purpose. In summary, the algorithm is summarized as follows, and the first step is to initialize the position. The flowchart is shown in Fig. 4.
The detailed pseudo code of the improved GWO algorithm is shown in algorithm 1.
Compared to the basic GWO algorithm, the improved grey wolf algorithm has the following characteristics: (1) Improved grey wolf algorithm introduces levy flight to balance global exploration and local mining capabilities; (2)Implement the mutation strategy for the current individual to enhance the diversity of the population and reduce the probability of the algorithm falling into local optimum.

Experimental results and analysis
The experimental platform is a notebook with 8G memory and Core i5. The program is implemented in PYTHON language and tested on public and self-made photo collections.
Through the evaluation of the camera model, the proposed improved grey wolf algorithm is implemented to obtain camera parameters. The 2-dimensional coordinates of the acquired image are as compared with the actual measured image coordinates to confirm the feasibility and superiority of the algorithm.
Step 1.To get the actual pixel coordinates in the image by Using PYTHON-OPENCV.
Step 2. Based on the parameters obtained by Zhang's method 7 , determine the upper and lower bounds of the iteration. In order to avoid slow convergence caused by too large optimization range, the obtained internal parameters ± 80 and distortion coefficient ± 10 are used as the upper and lower bounds of the optimization. Set the number of grey wolves in the wolf algorithm is 40, and it runs 400 times in total. The objective function is established, that is, the average value of the pixel distance between the actual coordinates and the optimized reprojection. Assuming we have N calibration points. Each point is identical in size and placed in the same noise environment. Objective function is established as follows: where N Denotes the number of calibrated corner points,p j represents the actual pixel coordinates of corner j , and p ′ j represents the reprojection coordinates of corner j . R and T are the rotation matric and translation vector corresponding to the jth image.
Step 3. Randomly generate an initial value between the upper and lower bounds, and the initial value is 9 dimensions.
Step 4. Calculate the value of f , obtain the wolves positions of α , β and γ according to the sorting result, and then optimize.
Step 5. Update the final position of the movement according to Eq. (14) and (15).
Step 6. If the iteration number reaches the maximum number, the iteration will be stopped; else return to step 4 to continue running.
Step 7. Finally, the best individual position and the best function value are returned, and the output is the mean value of the calibration error.
Experiments on public data sets. According to the article 11 , experiment in the public picture collection.
The series of photos are from the sample pictures of the calibration toolbox (http:// www. vision. calte ch. edu/ bougu etj/ calib_ doc/ htmls/ examp le. html). The size of the grid is 30mmx30mm. Each picture has 182 corner points. The unequal row-column dimension is adopted, which can always determine the direction of chessboard As can be seen from the Table 1, the average error value of the proposed method in this paper is 0.102 pixel, which is better than the 0.17 pixel in the paper 12 and improves the calibration accuracy to a certain extent.
In order to test the robustness of the algorithm, noise was added to the sample images and then compared with the calibration results of grey Wolf algorithm, particle swarm optimization algorithm and Zhang's calibration method. Examples and results of images with noise are shown in Figs. 5 and 6 respectively.
As can be observed in from Fig. 6, compared with the original image calibration results at the beginning, the average calibration error of the method in this paper has been increased by 0.07 pixels after adding noise, but the error is still the lowest in each method.
Two experiments show that the proposal algorithm improves the calibration accuracy effectively and has the characteristics of stability.     Fig. 7. Figure 8 shows the change of the objective function based on the hybrid algorithm in the process of 400 iteration. Table 2 lists the different results obtained by the four methods.
As can be seen from the figure above, the value of the objective function is in a stable state after about 260 calculations. Within 100 iterations, the algorithm converges rapidly and finds the extreme value in the local region from 100 to 250 times. It can be seen from Table 2 that the average error obtained by the four methods is all less than 0.1 pixel, which is due to the high accuracy of the experimental instrument in this paper. But relatively speaking, the calibration error value of this paper is 0.026 pixel, which is better than the 0.105 pixel   www.nature.com/scientificreports/ value of Zhang's calibration method and less than the 0.075 pixel value of grey wolf algorithm and the 0.058 pixel value of particle swarm optimization algorithm.

Conclusion and future works
Based on the idea of bionics, this paper establishes a gray wolf optimization algorithm model based on Levy flight and mutation mechanism. The algorithm is applied to the field of 3D reconstruction engineering. The nonlinear solution problem in the camera calibration process is solved. Through field experiments and data processing, the following conclusions are drawn: (1) Using the principle of bionics, the classical calibration algorithm is combined with the gray wolf algorithm. It effectively avoids the problem that the traditional calibration algorithm is easy to fall into the local optimal solution. The stability and robustness of the camera calibration algorithm are greatly improved. (2) Add Levy flight and mutation mechanism to the gray wolf optimization algorithm model. Optimize the convergence and population diversity of the algorithm. This can make the algorithm suitable for the actual camera calibration calculation, and effectively improve the calibration accuracy and stability. (3) Through the optimization algorithm model established in this paper, the parameter matrix of the camera is solved nonlinearly. The average error of reprojection of the calculation results is 0.026 pixels, which is better than that of PSO, GWO and other methods.
In the follow-up work, the research group will also conduct a detailed analysis of the reconstructed calibration point errors. The influence of the optimization algorithm model parameters on the calibration results is studied. The algorithm is then optimized and improved for a better accuracy and convergence speed of the camera calibration algorithm. In addition, the method in this paper has only been verified in the camera calibration calculation for the time being. In the future, we will continue to explore the feasibility of applying the algorithm in this paper to other fields.
It is worth mentioning that the work in this paper is helpful for the application of bionic optimization algorithms in the field of camera calibration, such as DMO, EOS, RSA, etc. mentioned above. This will provide a reference for future scholars' work to apply these excellent bionic algorithms to camera calibration.